{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Non-covalent interactions benchmark"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "start_time": "2023-05-01T23:12:06.304474Z",
     "end_time": "2023-05-01T23:12:06.663090Z"
    }
   },
   "outputs": [],
   "source": [
    "import pymolpro\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "start_time": "2023-05-01T23:12:06.664237Z",
     "end_time": "2023-05-01T23:12:06.665351Z"
    }
   },
   "outputs": [],
   "source": [
    "backend = 'local'  # If preferred, change this to one of the backends in your ~/.sjef/molpro/backends.xml that is ssh-accessible\n",
    "project_name = 'Non-covalent interactions benchmark'\n",
    "parallel = None  # how many jobs to run at once"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "start_time": "2023-05-01T23:12:06.666512Z",
     "end_time": "2023-05-01T23:12:06.667708Z"
    }
   },
   "outputs": [],
   "source": [
    "methods = {\"HF\": \"df-hf\", \"MP2\": \"df-mp2\", \"LMP2\": \"df-lmp2\", }\n",
    "bases = ['aug-cc-pVDZ', 'aug-cc-pVTZ', 'aug-cc-pVQZ']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "outputs": [],
   "source": [
    "db = pymolpro.database.load(\"GMTKN55_S22\").subset('small')"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-01T23:12:06.668910Z",
     "end_time": "2023-05-01T23:12:06.671134Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "outputs": [],
   "source": [
    "results = {}\n",
    "for method in methods:\n",
    "    results[method] = {}\n",
    "    for basis in bases:\n",
    "        results[method][basis] = pymolpro.database.run(db, methods[method], basis, location=project_name,\n",
    "                                                       backend=backend, parallel=parallel)"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-01T23:12:06.673321Z",
     "end_time": "2023-05-01T23:12:08.135297Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "outputs": [],
   "source": [
    "for method in methods:\n",
    "    for result in pymolpro.database.basis_extrapolate(results[method].values(),results['HF'].values()):\n",
    "        results[method][result.basis]=result\n",
    "    for basis in results[method]:\n",
    "        if basis not in bases: bases.append(basis)"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-01T23:12:08.136538Z",
     "end_time": "2023-05-01T23:12:08.138419Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "outputs": [
    {
     "data": {
      "text/plain": "             DF-HF      DF-MP2     DF-LMP2\n       aug-cc-pVQZ aug-cc-pVQZ aug-cc-pVQZ\nMSD          -6.08        0.22       -0.54\nSTDEVD        2.03        0.21        0.26\nMAD           6.08        0.22        0.54\nMAXD          7.52        0.42        0.79\nRMSD          6.30        0.28        0.58",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead tr th {\n        text-align: left;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr>\n      <th></th>\n      <th>DF-HF</th>\n      <th>DF-MP2</th>\n      <th>DF-LMP2</th>\n    </tr>\n    <tr>\n      <th></th>\n      <th>aug-cc-pVQZ</th>\n      <th>aug-cc-pVQZ</th>\n      <th>aug-cc-pVQZ</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>MSD</th>\n      <td>-6.08</td>\n      <td>0.22</td>\n      <td>-0.54</td>\n    </tr>\n    <tr>\n      <th>STDEVD</th>\n      <td>2.03</td>\n      <td>0.21</td>\n      <td>0.26</td>\n    </tr>\n    <tr>\n      <th>MAD</th>\n      <td>6.08</td>\n      <td>0.22</td>\n      <td>0.54</td>\n    </tr>\n    <tr>\n      <th>MAXD</th>\n      <td>7.52</td>\n      <td>0.42</td>\n      <td>0.79</td>\n    </tr>\n    <tr>\n      <th>RMSD</th>\n      <td>6.30</td>\n      <td>0.28</td>\n      <td>0.58</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 50,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.set_option('display.precision', 2)\n",
    "method_errors = pymolpro.database.analyse([results[method]['aug-cc-pVQZ'] for method in methods], db, 'kJ/mol')[\n",
    "                    'reaction statistics']\n",
    "with open(project_name + '.method_errors.tex', 'w') as tf:\n",
    "    tf.write('\\\\ifx\\\\toprule\\\\undefined\\\\def\\\\toprule{\\\\hline\\\\hline}\\n\\\\def\\\\midrule{\\\\hline}\\n\\\\def\\\\bottomrule{\\\\hline\\\\hline}\\\\fi') # or \\usepackage{booktabs}\n",
    "    tf.write(method_errors.style.format(precision=2).to_latex(hrules=True,multicol_align='c',caption='Method errors'))\n",
    "method_errors"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-02T12:26:47.371530Z",
     "end_time": "2023-05-02T12:26:47.382660Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "reaction energies\n",
      "        DF-HF      DF-MP2     DF-LMP2\n",
      "  aug-cc-pVQZ aug-cc-pVQZ aug-cc-pVQZ\n",
      "2    5.30e-03    8.11e-03    7.65e-03\n",
      "1    2.13e-03    5.09e-03    4.78e-03\n",
      "8   -5.94e-04    8.36e-04    7.38e-04\n",
      "reaction energy deviations\n",
      "        DF-HF      DF-MP2     DF-LMP2\n",
      "  aug-cc-pVQZ aug-cc-pVQZ aug-cc-pVQZ\n",
      "2   -2.65e-03    1.59e-04   -3.02e-04\n",
      "1   -2.87e-03    9.25e-05   -2.12e-04\n",
      "8   -1.43e-03   -3.85e-06   -1.02e-04\n",
      "reaction statistics\n",
      "             DF-HF      DF-MP2     DF-LMP2\n",
      "       aug-cc-pVQZ aug-cc-pVQZ aug-cc-pVQZ\n",
      "MSD      -2.32e-03    8.25e-05   -2.05e-04\n",
      "STDEVD    7.72e-04    8.18e-05    9.99e-05\n",
      "MAD       2.32e-03    8.51e-05    2.05e-04\n",
      "MAXD      2.87e-03    1.59e-04    3.02e-04\n",
      "RMSD      2.40e-03    1.06e-04    2.21e-04\n",
      "molecule energies\n",
      "          DF-HF      DF-MP2     DF-LMP2\n",
      "    aug-cc-pVQZ aug-cc-pVQZ aug-cc-pVQZ\n",
      "02      -152.14     -152.71     -152.71\n",
      "08       -80.43      -80.86      -80.85\n",
      "01      -112.45     -112.96     -112.96\n",
      "01a      -56.22      -56.48      -56.48\n",
      "02a      -76.07      -76.35      -76.35\n",
      "08a      -40.22      -40.43      -40.43\n"
     ]
    }
   ],
   "source": [
    "for key, table in pymolpro.database.analyse([results[method]['aug-cc-pVQZ'] for method in methods], db).items():\n",
    "    print(key)\n",
    "    print(table)"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-02T06:43:48.649242Z",
     "end_time": "2023-05-02T06:43:48.652936Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "outputs": [
    {
     "data": {
      "text/plain": "        DF-HF      DF-MP2     DF-LMP2\n  aug-cc-pVQZ aug-cc-pVQZ aug-cc-pVQZ\n2       -6.96        0.42       -0.79\n1       -7.52        0.24       -0.56\n8       -3.76       -0.01       -0.27",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead tr th {\n        text-align: left;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr>\n      <th></th>\n      <th>DF-HF</th>\n      <th>DF-MP2</th>\n      <th>DF-LMP2</th>\n    </tr>\n    <tr>\n      <th></th>\n      <th>aug-cc-pVQZ</th>\n      <th>aug-cc-pVQZ</th>\n      <th>aug-cc-pVQZ</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>2</th>\n      <td>-6.96</td>\n      <td>0.42</td>\n      <td>-0.79</td>\n    </tr>\n    <tr>\n      <th>1</th>\n      <td>-7.52</td>\n      <td>0.24</td>\n      <td>-0.56</td>\n    </tr>\n    <tr>\n      <th>8</th>\n      <td>-3.76</td>\n      <td>-0.01</td>\n      <td>-0.27</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.set_option('display.precision', 2)\n",
    "method_errors = pymolpro.database.analyse([results[method]['aug-cc-pVQZ'] for method in methods], db, 'kJ/mol')[\n",
    "                    'reaction energy deviations']\n",
    "method_errors"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-02T12:25:56.377711Z",
     "end_time": "2023-05-02T12:25:56.388995Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "outputs": [
    {
     "data": {
      "text/plain": "           DF-LMP2                                                      \n       aug-cc-pVDZ aug-cc-pVTZ aug-cc-pVQZ aug-cc-pV[23]Z aug-cc-pV[34]Z\nMSD          -1.21       -4.05       -0.54          -4.85           2.07\nSTDEVD        1.36        5.71        0.26           7.71           3.80\nMAD           1.44        4.05        0.54           4.85           2.24\nMAXD          2.23       10.62        0.79          13.76           6.46\nRMSD          1.64        6.17        0.58           7.95           3.73",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead tr th {\n        text-align: left;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr>\n      <th></th>\n      <th colspan=\"5\" halign=\"left\">DF-LMP2</th>\n    </tr>\n    <tr>\n      <th></th>\n      <th>aug-cc-pVDZ</th>\n      <th>aug-cc-pVTZ</th>\n      <th>aug-cc-pVQZ</th>\n      <th>aug-cc-pV[23]Z</th>\n      <th>aug-cc-pV[34]Z</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>MSD</th>\n      <td>-1.21</td>\n      <td>-4.05</td>\n      <td>-0.54</td>\n      <td>-4.85</td>\n      <td>2.07</td>\n    </tr>\n    <tr>\n      <th>STDEVD</th>\n      <td>1.36</td>\n      <td>5.71</td>\n      <td>0.26</td>\n      <td>7.71</td>\n      <td>3.80</td>\n    </tr>\n    <tr>\n      <th>MAD</th>\n      <td>1.44</td>\n      <td>4.05</td>\n      <td>0.54</td>\n      <td>4.85</td>\n      <td>2.24</td>\n    </tr>\n    <tr>\n      <th>MAXD</th>\n      <td>2.23</td>\n      <td>10.62</td>\n      <td>0.79</td>\n      <td>13.76</td>\n      <td>6.46</td>\n    </tr>\n    <tr>\n      <th>RMSD</th>\n      <td>1.64</td>\n      <td>6.17</td>\n      <td>0.58</td>\n      <td>7.95</td>\n      <td>3.73</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.set_option('display.precision', 2)\n",
    "basis_errors = pymolpro.database.analyse([results['LMP2'][basis] for basis in bases], db, 'kJ/mol')[\n",
    "    'reaction statistics']\n",
    "with open(project_name + '.basis_errors.tex', 'w') as tf:\n",
    "    tf.write('\\\\ifx\\\\toprule\\\\undefined\\\\def\\\\toprule{\\\\hline\\\\hline}\\n\\\\def\\\\midrule{\\\\hline}\\n\\\\def\\\\bottomrule{\\\\hline\\\\hline}\\\\fi') # or \\usepackage{booktabs}\n",
    "    tf.write(basis_errors.style.format(precision=2).to_latex(hrules=True,multicol_align='c',caption='Basis errors'))\n",
    "basis_errors"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-02T12:23:16.813556Z",
     "end_time": "2023-05-02T12:23:16.827929Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "outputs": [
    {
     "data": {
      "text/plain": "<Figure size 1800x600 with 3 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABc4AAAIlCAYAAADsVGvuAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABe2UlEQVR4nO3deZRcdZk//qeql+rukI3sgZAEXFjCJjgs40iCGkCBYREIOkIQUEAWQX44GReWESOyyBdEFgeC7KjIYVWEYREFZBFEECOMwUAWQtbO2lvV7w9M20227qS7qu6t1+ucOid161bdp6rS/XS976c+n0yhUCgEAAAAAAAQERHZUhcAAAAAAADlRHAOAAAAAAAdCM4BAAAAAKADwTkAAAAAAHQgOAcAAAAAgA4E5wAAAAAA0IHgHAAAAAAAOhCcAwAAAABAB4JzAAAAAADoQHAOJN4999wTmUwmrrnmmnXu8/DDD0cmk4nLLrssMplMly4R0aX9zjvvvCI9UwBIpu726tUKhULcdtttse+++8bAgQOjrq4uttlmmzjttNNi1qxZne7/+OOPd6lvP/744731NAGgLG1MHx4zZkx77zz11FPb91u+fHlMmjQpPvzhD0ffvn2jT58+scMOO8R3vvOdWL58+Xrr+OY3vxmZTCbGjRu3xm277LJL+/EOPPDA9u0d61jXZfLkyd18RaBrMoVCoVDqIgA2RWtra4waNSpGjRoVzz777Fr3+dznPhc///nP4/XXX485c+Z0uu3QQw+NbbbZJi655JJO2/fcc8945pln1nnMY445JmbNmhVPPvlk/Mu//EvPPBkASKHu9OpZs2bFkCFDIp/Px+c+97m488474+ijj46jjjoq+vfvHy+//HJcfPHFsXLlynjwwQfbe3BjY2P8+c9/Xutjv/322/H5z38+tthii3jhhRdi4MCBvfZcAaDcbEwfHjNmTOywww7xrW99K4YPHx5jxoyJiIjFixfHl7/85dh3331j7Nixkc1m4ze/+U1873vfi49//OPxyCOPrPXxX3rppdhzzz1jwIABMXjw4HjllVc63f6nP/0pli9fHoceemjstttucf/990dExIsvvhhNTU1rfczTTz89nn/++fjFL34RhxxyyMa9OLAe1aUuAGBTVVdXxzHHHBPf//7345VXXlnj7PXixYvj7rvvjoMPPjhGjx4do0eP7nR7LpeLAQMGxJ577rnGY69tW8R7DXrGjBlx7bXXCs0BYAO606uHDBkSEREXXXRR3HnnnfG9730vvv71r7fvO378+DjqqKNijz32iMMOOyxee+216Nu3b/Tr12+tfbupqSnOOOOMqKqqil/84hdCcwAqzsb04YiIIUOGrNFbBwwYEHfeeWenbZ/85Cejqakpvv/978ff/va32HrrrTvd3traGscdd1x8+ctfjj/+8Y8xf/78NWrccccdI+K9z+cd7brrrmt9Tpdddlk899xzMWXKFKE5vcZULZAQb7zxRhx33HHxwQ9+MBoaGmKLLbaIgw46KP70pz912u/GG2+MTCYTb775Zqftq7++3PHryYVCIb773e/G6NGjo66uLnbfffd4+OGHY/z48TF+/Pgu1dXU1BQXXHBBbLfddlFXVxeDBg2KCRMmxFNPPdW+Tz6fjyuvvDJ22WWXqK+vbw+p77333g0+/vjx42PcuHHx5JNPxp577hn19fWxxRZbxLe+9a1oa2tr3+/444+PiIhp06at8Ri33357rFq1Kr74xS926TltyM033xxXXnllHH/88fGlL32pRx4TgOTTq3uuVzc3N8fFF18c2223XZxzzjlr7D9s2LCYOnVqzJo1a62P19Epp5wSzz77bFx33XWxyy67bPD5AJBM+nBpPzOvDtyrq9cco/u9730vFi5cGBdeeOFGP35Hjz32WHz961+PiRMnxne+850eeUxYG8E5JMTs2bNj0KBB8b3vfS9+9atfxVVXXRXV1dWxxx57xPTp0zfqMb/xjW/EN77xjdh///3jnnvuiZNOOilOOOGE+Otf/9ql+7e2tsYBBxwQ//3f/x0HHnhg3H333XHjjTfG3nvvHTNnzmzfb/LkyXHGGWfERz/60bjzzjvjjjvuiIMPPniNP1TWZe7cuTFp0qT4/Oc/H/fcc0989rOfje985ztxxhlntO/zoQ99KD72sY/FLbfcEi0tLZ3uP23atNhiiy1iv/3269Lx1ufFF1+ML3/5y/HRj340rrrqqk1+PADSQ6/uuV79wgsvxKJFi+Lggw9uX3fk/Q466KDIZrPx0EMPrbOuq6++Om644YY47bTT4j/+4z+69FwASCZ9uLifmQuFQrS2tkZjY2P86le/iksvvTSOPvro2GqrrTrt9+c//zm+853vxNVXXx2bbbZZlx9/XWbOnBlHHXVUbLnllnH77bdHNivapBcVgERqbW0tNDc3Fz74wQ8WzjzzzPbt06ZNK0REYcaMGZ32f+yxxwoRUXjssccKhUKhsHDhwkIulyscddRRnfZ7+umnCxFR2GeffTZYw0033VSIiMKPf/zjde7zm9/8phARhW984xtdfm4d7bPPPoWIKNxzzz2dtp944omFbDZb+Pvf/96+bfVz/8UvftG+7ZVXXtng8UePHl34zGc+s8Fa3n333cLo0aMLQ4YMKcycOXMjng0AlUSv3vhefccddxQionDNNdes99jDhg0r7LDDDmu97Xe/+12hpqam8G//9m+F5ubmjXlqACSYPrxpn5lHjx5dOPbYY9d53Ntvv70QEe2X4447rtDS0tJpn7a2tsIee+xROProozvVu67evfq46/p8vnLlysJuu+1WqK+vL/zhD39Y52NAT3FaBhKitbU1vvvd78b2228ftbW1UV1dHbW1tfH666/Ha6+91u3He+aZZ6KpqSmOPPLITtv33HPP9kU/Vmtra4vW1tb2Sz6fj4iIX/7yl1FXV7fer3P98pe/jIiIr3zlK+vcJ5/Pd3r8jl8ni4jo27dvHHzwwZ22fe5zn4t8Ph+/+c1v2rcdeeSR0bdv37jhhhvat91www2RyWTiuOOOW+fxu6KtrS0mTZoUb7/9dtx5550xatSoTXo8ANJHry5+ry4UCmsdkT5nzpz47Gc/G0OGDImf/vSnUVNT0+3HBiBZ9OHi9uH99tsvnnvuuXj00UfjwgsvjLvuuisOP/zw9uce8d485K+//npcfvnl3XrsdTnppJPihRdeiGuvvXadc59DTxKcQ0KcddZZ8a1vfSsOOeSQuO++++L3v/99PPfcc7HzzjvHypUru/14CxYsiIj35gh9v/dv22abbaKmpqb9csEFF0RExLvvvhsjR45c71ej3n333aiqqorhw4evc58LLrig0+Nvs802660nItofb/XziIhoaGiISZMmxa9+9auYO3dutLa2xi233BL77LPPGo/ZXeecc0787//+b1x00UUxYcKETXosANJJr+5sU3r16q95z5gxY501LV++PObPn7/Gyezm5uY4/PDDY8GCBfHzn/98vc8LgPTQhzvr7c/MAwcOjN133z0mTJgQ//Vf/xXXXXdd3HvvvXHPPfdExHtTqnz729+Oc889N2pra2Px4sWxePHi9hMLixcv7tb7cuWVV8ZPfvKTOPXUU+MLX/hCt2qFjbXmjP1AWbrlllvimGOOie9+97udts+fPz8GDBjQfr2uri4i3luA5P37dTRo0KCIiHjnnXfWONbcuXM7nUG/7777Oj3eyJEjI+K9xT9++9vfRj6fX+cfAkOGDIm2traYO3dujBgxYq37fOlLX4oDDzyw/fr7V9FeV40dn8dqxx9/fPz4xz+Om266KT70oQ/FvHnz4tJLL13rcbvq9ttvj8suuyyOOuqo+NrXvrZJjwVAeunVa9bY8Xms1pVevdtuu8Xmm28e9957b0ydOnWto8rvvffeyOfzse+++3baftppp8XTTz8dP/rRj2KvvfZa6/MBIH304TVr7Pg8VuuNz8wREf/yL/8SEdE+//vf/va3WLlyZZxxxhmd5lpfbeDAgXHGGWd0aTT6k08+GV/72tfiYx/7WFx22WWbXCt0lRHnkBCZTGaN5vjAAw/ErFmzOm1b3bxffvnlTtvfvxr3HnvsEblcLu68885O25955pn4+9//3mnbjjvuGLvvvnv7ZfUfAQcccECsWrUqbrzxxnXWfcABB0TEe4tzrcvIkSM7Pf6OO+7Y6falS5euUf9tt90W2Ww2Pv7xj6/xvMaNGxfTpk2LadOmRf/+/ePwww9f57E35OWXX44TTjghxo0bF9dff/1GPw4A6adX91yvrq2tjf/v//v/4rXXXouLL754jXrmzZsXU6ZMiQEDBsTkyZPbt//P//xPXHfddXHcccfFySefvM7nA0D66MOl+cy82mOPPRYRER/4wAciImKXXXaJxx57bI3LzjvvHGPGjInHHnssTj311A0+7qxZs+KII46IwYMHx89+9jPTr1FURpxDQhx44IFx4403xrbbbhs77bRTvPDCC3HxxRfHlltu2Wm/j370o/HhD384zj777GhtbY2BAwfG3XffHb/97W877bf55pvHWWedFVOnTo2BAwfGoYceGm+//Xacf/75MWLEiC6tTH300UfHtGnT4qSTTorp06fHhAkTIp/Px+9///vYbrvtYtKkSfFv//Zv8YUvfCG+853vxDvvvBMHHnhg5HK5ePHFF6OhoSFOO+20DR5n0KBBcfLJJ8fMmTPjQx/6UDz44IPx4x//OE4++eQ1VuyOiPjiF78YZ511VkyfPj2+/OUvR319/QaPsTaLFi2KQw45JJqamuLrX/96/OlPf1rrfkOGDNnkqWAASD69umd79TnnnBMvvfRSfP3rX48//vGPcdRRR0X//v3j5ZdfjosvvjjeeeeduP/++2Pw4MEREfHss8/GqaeeGsOHD49jjjkmnnnmmbXWus0228SQIUM2+JwASBZ9uDifma+99tp48sknY+LEiTFq1KhYvnx5PPnkk3HllVfG3nvvHf/+7/8eEREDBgyI8ePHr3H/AQMGRGtr61pve7/m5uY47LDD4p133olLL7003nzzzXjzzTfX2K9fv36x/fbbd/k5QJeVenVSoGsWLVpUOP744wtDhw4tNDQ0FD72sY8VnnzyycI+++yzxmref/3rXwsTJ04s9OvXrzBkyJDCaaedVnjggQc6rRBeKBQK+Xy+8J3vfKew5ZZbFmpraws77bRT4f777y/svPPOhUMPPbRLda1cubLw7W9/u/DBD36wUFtbWxg0aFBh3333LTz11FPt+7S1tRV+8IMfFMaNG1eora0t9O/fv7DXXnsV7rvvvg0+/uoVtx9//PHC7rvvXsjlcoURI0YU/uu//muNFbtXe/fddwu1tbWFiCg8++yzGzzGulbtXr2q+oYu61tpHIDKoVf3fK/O5/OFm2++ubDPPvsU+vfv3957P/zhDxdee+21Tvuee+65Xerb06ZN69LrBkCy6MM924dHjx691s+6v/vd7woHHnhgYeTIkYXa2tpCQ0NDYeeddy7893//d2H58uVdrnddOn4+nzFjRpd6+/vfX+gpmUKhUOjlbB5IkBkzZsS2224b5557bvzXf/1XqcuJ8ePHx/z58+OVV14pdSkAUBYqvVefcMIJ8ZOf/CTuuuuuOPjgg4tyTABYrVL68JgxY2KfffaJ66+/PrLZbJdG2G+Ktra2KBQK8YEPfCDGjRsX999/f68eD7rCVC1Qwf74xz/G7bffHnvvvXf069cvpk+fHt///vejX79+cfzxx5e6PACoeHr1mq699tp455134sgjj4z77rsvPvWpT5W6JABSqtL78E033RQ33XRTfOUrX4kf/vCHvXqs3XbbLf74xz9GRMS4ceN69VjQVYJzqGB9+vSJ559/Pq6//vpYvHhx9O/fP8aPHx8XXnhhDBs2rNTlAUDF06vXVFVVFffdd1+pywCgAlRyH77vvvuiqakpIiKGDh3a68e77bbbYsWKFRHx3jzoUA5M1QIAAAAAAB307gRFAAAAAACQMIJzAAAAAADoQHAOAAAAAAAdWBy0B+Tz+Zg9e3b07ds3MplMqcsBIGUKhUIsXbo0Ro4cGdmsc94bS78GoLfp2T1DzwagN3W1XwvOe8Ds2bNj1KhRpS4DgJR76623Yssttyx1GYmlXwNQLHr2ptGzASiGDfVrwXkP6Nu3b0S892L369evxNUAkDaNjY0xatSo9n7DxtGvAehtenbP0LMB6E1d7deC8x6w+qtj/fr109QB6DW+qrxp9GsAikXP3jR6NgDFsKF+bdI1AAAAAADoQHAOAAAAAAAdCM4BAAAAAKADwTkAAAAAAHQgOAcAAAAAgA4E5wAAAAAA0IHgHAAAAAAAOhCcAwAAAABAB4JzAAAAAADoQHAOAAAAAAAdCM4BAAAAAKADwTkAAAAAAHQgOAcAAAAAgA4E5wAAAAAA0IHgHAAAAAAAOhCcAwAAAABAB9WlLgAAAAAAIElWNLeWuoSSaKitnDi5cp4pAAAAAEAP2P7bD5W6hJJ483ufKXUJRSM4BwAAAABYiyUrW+LVWUtKXUbZeOqN+Wts23rIZjG8f10JquldgnMAAAAAgLVoXNkSK5rb1th+1dG7Fr2Wpta2OOtnL0dExGVH7BS56qqi17C212LxymbBOQAAAABApVjWtPa5zHM1xQ+tOx2/uqrkNay2bFU653vPlroAAAAAAIBytK7gnH9a2tQahUKh1GX0OME5AAAAAMD7FAqF1I6m7kltbYVY1ZIvdRk9TnAOAAAAAPA+K1vaoi2fvpHUvWFpU0upS+hxgnMAAAAAgPcx2rzrlqbwtRKcAwAAAAC8z1Lzm3dZGk8yCM4BAAAAAN4njWFwb1mewpMMgnMAAAAAgPdZlsIwuLesaG6L1rZ0LRAqOAcAAAAA6KC1LR8rm9tKXUaiLG9K1+slOAcAAAAA6MBo8+5b2tRS6hJ6lOAcAAAAAKADwXn3GXEOAAAAAJBigvPuW2bEOQAAAABAei1bJTjvrqUpe80E5wAAAAAAHRhx3n2tbYVY1ZKe6VoE5wAAAAAA/7CqpS1a2wqlLiOR0nTCQXAOAAAAAPAPaQp/i215il47wTkAAAAAwD+Y33zjpWmec8E5AAAAAMA/GHG+8dL02gnOAQAAAAD+IU3hb7GtaG6NQiEd88MLzgEAAAAAIiKfL8SKZsH5xsrnI1Y0t5W6jB4hOAcAAAAAiIgVLW2Rz5e6imRLy4h9wTkAAAAAQEQsT0noW0ppWSBUcA4AAAAAEOkJfUspLScfBOcAAAAAAJGe0LeUTNUCAAAAAJAiaQl9S2llc1u0tiV/onjBOQAAAABQ8Vrb8rGyua3UZaTC8qbkv46CcwAAAACg4qUh7C0Xy5qTP3JfcA4AAAAAVLylTS2lLiE1lqVgkVXBOQAAAABQ8Yw47zlpmCtecA4AAAAAVLxlRpz3GMF5mTnvvPMik8l0ugwfPny993niiSdit912i7q6uth6663jmmuuKVK1AAAAAEC5WJqC6UXKRUtrPla1JHsEf3WpC+hpO+ywQzzyyCPt16uqqta574wZM+LTn/50nHjiiXHLLbfE7373uzjllFNiyJAhcfjhhxejXAAAAACgxFa1tEVrW6HUZaTK8qbWqKtZdzZb7lIXnFdXV29wlPlq11xzTWy11VZx+eWXR0TEdtttF88//3xccsklgnMAAAAAqBBpmFqk3Cxrao1Bm+VKXcZGS9VULRERr7/+eowcOTLGjh0bkyZNir/97W/r3Pfpp5+OiRMndtq23377xfPPPx8tLeue06ipqSkaGxs7XQCA8qJfA0Ay6NlAOVhmmpYel/Spb1IVnO+xxx5x0003xUMPPRQ//vGPY+7cubH33nvHggUL1rr/3LlzY9iwYZ22DRs2LFpbW2P+/PnrPM7UqVOjf//+7ZdRo0b16PMAADadfg0AyaBnA+XAiPOel/TXNFXB+QEHHBCHH3547LjjjvHJT34yHnjggYiI+MlPfrLO+2QymU7XC4XCWrd3NGXKlFiyZEn75a233uqB6gGAnqRfA0Ay6NlAOUj66OhytKK5NfL55M4bn7o5zjvq06dP7LjjjvH666+v9fbhw4fH3LlzO22bN29eVFdXx6BBg9b5uLlcLnK55M7PAwCVQL8GgGTQs4FSy+cLsaJZcN7T8vmI5c2t0beuptSlbJRUjTh/v6ampnjttddixIgRa719r732iocffrjTtl//+tex++67R01NMt9QAAAAAKDrljW3RiG5A6PLWpKna0lVcH722WfHE088ETNmzIjf//738dnPfjYaGxvj2GOPjYj3vv51zDHHtO9/0kknxd///vc466yz4rXXXosbbrghrr/++jj77LNL9RQAAAAAgCKyMGjvSfJrm6qpWt5+++04+uijY/78+TFkyJDYc88945lnnonRo0dHRMScOXNi5syZ7fuPHTs2HnzwwTjzzDPjqquuipEjR8YVV1wRhx9+eKmeAgAAAABQREkeFV3ulib4tU1VcH7HHXes9/Ybb7xxjW377LNP/OEPf+iligAAAACAcmZh0N6T5Nc2VVO1AAAAAAB0x9JVLaUuIbVaWvOxqqWt1GVsFME5AAAAAFCRVrW0RWublUF7U1KnwhGcAwAAAAAVKclTiSRFUhcIFZwDAAAAABXJNC29L6knJwTnAAAAAEBFSmqomyRJPTkhOAcAAAAAKlJS599OkhXNbdHali91Gd0mOAcAAAAAKk5LWz5WNreVuoyKkMQTFIJzAAAAAKDimKaleJL4WgvOAQAAAICKk9S5t5OoMYGvteAcAAAAAKg4SRwFnVRJfK0F5wAAAABAxUniKOikWt7UGvl8odRldIvgHAAAAACoKG35QqxosjBosRQKEUsTtkCo4BwAAAAAqCjmNy++pL3mgnMAAAAAoKI0rkzW6Oc0SNprLjgHAAAAACqK+c2Lz4hzAAAAAIAytnRVskY/p8Hy5mQtECo4BwAAAAAqRmtbPpYnbKHKNMjnk7VAqOAcAAAAAKgYRpuXTuPK5EzXIjgHAAAAACqG4Lx0kjS3vOAcAAAAAKgYSQpv0yZJJy0E5wAAAABAxUjSdCFps7ypNdoSskCo4BwAAAAAqAgtbflY0dxW6jIqVqEQsTQhI/4F5wAAAABARTDavPQaVyZjuhbBOQAAAABQERoTNMd2Wi1JyMkLwTkAAAAAUBGSEtqmWVIWZxWcAwAAAAAVwVQtpbeyuS2aW/OlLmODBOcAAAAAQOqtaklGYFsJkjDyX3AOAAAAAKReEsLaSpGE6VoE5wAAAABA6gnOy0cS3gvBOQAAAACQeuY3Lx+NK1uiUCiUuoz1EpwDAAAAAKmWzxcSMT1IpWhtK8SK5rZSl7FegnMAAAAAINWWNbdG3rqgZaXcp2sRnAMAAAAAqbZkRXmHtJWo3L8BIDgHAAAAAFKt3Ec3V6JyP5khOAcAAAAAUs3CoOVnWVNrtOXLd4FQwTkAAAAAkFrNrfmyX4iyEhUK5X1CQ3AOAAAAAKRWuc+lXcnKeQodwTkAAAAAkFqLy3wu7UomOAcAAAAAKAEjzsuX4BwAAAAAoMgKhUJZh7OVrrk1HyvLdP55wTkAAAAAkErLmlqjra1Q6jJYj3I9sSE4BwAAAABSqVxDWf6pXN8jwTkAAAAAkErlGsryT+X6HgnOAQAAAIBUKtdQln9auqol2vLlN52O4BwAAAAASJ2WtnysaCrPhSf5p0LhvfC83AjOAQAAAIDUMdo8OcrxvRKcAwAAAACpU45hLGtXju+V4BwAAAAASJ1yDGNZu3J8rwTnAAAAAECqFAqFsgxjWbumlnysaimv+egF5wAAAABAqqxobou2tkKpy6Abyu1ER6qC86lTp8ZHP/rR6Nu3bwwdOjQOOeSQmD59+nrv8/jjj0cmk1nj8pe//KVIVQMAAAAAPancQlg2rNzes1QF50888UR85StfiWeeeSYefvjhaG1tjYkTJ8by5cs3eN/p06fHnDlz2i8f/OAHi1AxAAAAANDTyi2EZcPK7T2rLnUBPelXv/pVp+vTpk2LoUOHxgsvvBAf//jH13vfoUOHxoABA3qxOgAAAACgGMothGXDlq5qiXy+ENlsptSlRETKRpy/35IlSyIiYvPNN9/gvrvuumuMGDEiPvGJT8Rjjz3W26UBAAAAAL2gtS0fy5taS10G3ZTPRywto/ctVSPOOyoUCnHWWWfFxz72sRg3btw69xsxYkRcd911sdtuu0VTU1PcfPPN8YlPfCIef/zxdY5Sb2pqiqampvbrjY2NPV4/ALBp9GsASAY9G+hpS1e1RsG6oInUuLIl+tfXlLqMiEhxcH7qqafGyy+/HL/97W/Xu9+HP/zh+PCHP9x+fa+99oq33norLrnkknUG51OnTo3zzz+/R+sFAHqWfg0AyaBnAz3NNC3JtWRlS4wqdRH/kMqpWk477bS4995747HHHostt9yy2/ffc8894/XXX1/n7VOmTIklS5a0X956661NKRcA6AX6NQAkg54N9LTGVYLzpGoso5MeqRpxXigU4rTTTou77747Hn/88Rg7duxGPc6LL74YI0aMWOftuVwucrncxpYJABSBfg0AyaBnAz3NiPPkWtHcFi1t+aipKv1471QF51/5ylfitttui3vuuSf69u0bc+fOjYiI/v37R319fUS8dyZ71qxZcdNNN0VExOWXXx5jxoyJHXbYIZqbm+OWW26Ju+66K+66666SPQ8AAAAAoPtWtbRFU0u+1GWwCRpXtsSgzUp/QjVVwfnVV18dERHjx4/vtH3atGkxefLkiIiYM2dOzJw5s/225ubmOPvss2PWrFlRX18fO+ywQzzwwAPx6U9/ulhlAwAAAAA9wDQtyde4qlVw3tMKXVgu98Ybb+x0/ZxzzolzzjmnlyoCAAAAAIqlcWVrqUtgE5XLPOelnywGAAAAAKAHGHGefOXyHgrOAQAAAIBUWLrKiPOka2rJR1NrW6nLSNdULQAAAABAZVrV0hYtrcVZGLSppfjBbscwuVTBcq6mqijHWbqqNXKbFedY6yI4BwAAAAASr5hzY3/l9heLdqy1OetnL5fkuP9zzO5FOU7jypYYXOIFQk3VAgAAAAAkXqNpWlKjHKbcMeIcAAAAAEi8pUVcVPKqo3ct2rFWa2ptax9pftkRO0WuurRTmfQmwTkAAAAAQA8oZtharLm+13n86qqS19CbVrW0RXNrPmqrSzdhiqlaAAAAAIBEa2p9L2glPZY1lXbUueAcAAAAAEi0ZWUwtQc9q9TvqeAcAAAAAEi0Uo9Opuc1FnHO+rURnAMAAAAAiVYOi0nSs0p9MkRwDgAAAAAkWqlDVnreiubWKBQKJTu+4BwAAAAASKx8vhArmgXnaZPPR6xobivZ8QXnAAAAAEBirWxpi3y+1FXQG5aX8JsEgnMAAAAAILFKGa7Su0o5BY/gHAAAAABILPObp9fyJlO1AAAAAAB0WynnwaZ3LS/h3PWCcwAAAAAgsUzVkl4rmlujUCiU5NiCcwAAAAAgsYw4T698PqKptTQrvwrOAQAAAIBEWtXSFm350oxIpjhK9Y0CwTkAAAAAkEirWow2T7tSfaNAcA4AAAAAJJJpWtJvZYlOjgjOAQAAAIBEEpyn30ojzgEAAAAAus5ULelnxDkAAAAAQDcIztOvVO+x4BwAAAAASKRVLflSl0Ava20rREtb8d9nwTkAAAAAkDiFQiGaWo04rwRNrYJzAAAAAIANamrNR6FQ6ioohlJM1yI4BwAAAAASpxSjkCkNI84BAAAAALrANC2Vo1lwDgAAAACwYaUIUykNwTkAAAAAQBcIzitHKd7r6qIfkbK2orm11CWUREOtHwUAAACAJGlpszJopWhuK/60PNJCOtn+2w+VuoSSePN7nyl1CQAAAAB0Q0ubEeeVohQnSQTnKZfPF+LdZU0xe/HKWLSiudTllK1H//LOem/PVVfF8P51MbJ/fdTXVhWpKgAAAADWpVlwXjFKcZJEcJ5SS1e1xOzFq2LOkpXR2o0zMlcdvWsvVrV2Ta1tcdbPXo6IiMuO2Cly1cUPpvMb+Nlb2dwWM95dHjPeXR4D+9TGFgPqY2jfXGSzmeIUCAAAAEAn3cm8SLZSvNeC8xRpyxdi3tJV8failbFkRctGPUauprSjqXPVVSWvYUMWLW+ORcubo6Y6GyP718UWA+vNkQ4AAABQZK0bGglJarTlC1EoFCKTKd4gVmlfCqxobo23F62M2Yu7N7qcTdPSmo+/L1gRf1+wIjbfrDa2HFgfQzbLFfUHGAAAAKBSteXlYJWkNV+ImirBORtQKLw3d/nbi1bGwmXmLi+1hcuaY+Gy5sjVZGOLAfUxckB91JX5yHkAAACAJGsVnFeUtnwhihm3Cc4Tpqm1LWYvXhVvL1oRTS2+jlJumlry8bd3l8eM+ctjaN+62HJgfQzsU1vqsgAAAABSJy84ryjF/oaB4DwhlqxsibcWroh5S1dtcCHLTdHU0tZ7D76uY7a2rfXfxdTT86oXChHvNK6KdxpXxWZ11TFq84YY3q8uqiwmCgAAALDJCoVCFOTmFSVf5DdccF7mFixrijcXLI9Fyzdusc/u+srtLxblOOty1s9eLslx/+eY3XvtsZetao3XZjfGG/OWxVabN8SWA+ujpirba8cDAAAASDuDzStPsdeCFZyXqcZVLfH6O0uLFpjT+1pa8/F/85bF3xcsj60HbxajNq+3kCgAAADARij26GNKrxBGnFe8mQtWxOvzlpbk6yZXHb1r0Y/Z1NrWPtL8siN2ilx1uhfVbG0rxF/fWRrvLF0VO285IGqrjT4HoHtWNLeWuoSSaKj1pysAAO+Rm1eeYr/nPn2UmQXLmuKv7ywt2fF7eq7vbh+/uqrkNRTLkhUt8ersJbHrVgNLXQoAvSCfL8Ty5tZobev5v+52/e+He/wxk+DFb32qxx8zm83EZrlq65AAAECZK/a5EsF5mWk1QVNFKfZqwAD0npa2fCxa0RyLV7TEkpUtsXRVS9Hn4Eu7F/6+qFceN5OJ6JOrjgENNTGgvjYG9qlJ/TfgAAAgaYo91EVwXmaG9auLZUNa4835y0vylZOmlrbiH7O1ba3/LqZSjHLvW1cdO4zsX/TjAtBzWtvyMWfJqpi3dFUsXtFStN5tarWeVSi8t5j3slWt8XasjIj3+vSwfnUxYkBdqp4rAADQNYLzMrTNkM1iRP+6+PuCFTG3cVW09cJXvNflK7e/WLRjrc3qD+TF9j/H7F60Yw1oqIlRmzfE0L45i4MCJNicJStj+tylvTIVy4aUelqzSphabemq1li6aln8bf6yGDt4sxg7uE+pSwIAoAORSuUp9nsuOC9TDbXVsd2IfvGhYX1jwfKmaGrxXe/e9OHhfXv9GJlMxOZ9ai1sBpACrW35eHVWY8mO7xtixZPPR/zfvGUxeLPa6FtXU/TjAwCwdlnJecUp9gBUCV6Zq8pmYmjfuqId788X7Fe0Y5UTYTYA3VFdlY0tBtbHrEUrS3J83xArrkGb1UYffysA0E0rmltLXULJ+IxNMVjbvfIU+z33m4xONDcA6JrtRvSLQX1q4415y2JFc2lGYNO7cjXZGDOoT2w5sN70agAJVSgU4t2lTdHcVvxvcX/soseKfsxy8duvTyj6MetrqmLQZrmiH5fSyWQykc2+9w1BKkNVkZPzVKakP/rRj+Liiy+OOXPmxA477BCXX355/Nu//ds693/iiSfirLPOildffTVGjhwZ55xzTpx00klFrBgASKKh/epiSN9czF/WHLMWr4wFy5qKskCoxUF718A+NTFyQH0M61sXWUOZABLttTlLY/bi0nxDrJL9Zc7Skhz3A0M3izHWJako2Uwm8lH8NYcojWJPz5O64PzOO++Mr371q/GjH/0o/vVf/zWuvfbaOOCAA+LPf/5zbLXVVmvsP2PGjPj0pz8dJ554Ytxyyy3xu9/9Lk455ZQYMmRIHH744SV4BgBAkmQymRjSNxdD+uaiuTUf85auinlLm2LxiuZeG/1S6oU507Y4aCYT0b++Job2rYuh/XJRl6LnBlDpRg9qiKbWtliwrLnox77siJ2Kfsym1nxMufuViIiYeui4yFVni15DKWSzEUP71sWIAcWb6pbyUJXNRGub4LxSGHG+iS677LI4/vjj44QTToiIiMsvvzweeuihuPrqq2Pq1Klr7H/NNdfEVlttFZdffnlERGy33Xbx/PPPxyWXXCI4BwC6pbY6G1sObIgtBzZEa1s+Fq5ojoXL37usaOq56VwsDrrp6mqqYmCfmhjUJxeb96mN2goJFgAqTZ9cdey61cBY1dIW8xqbYv7y3j253VGp1gRZbXWAXgrFWJekqioTmzfUxuC+uRiyWU4vr1DV2Ww0hblaKkW14HzjNTc3xwsvvBD/+Z//2Wn7xIkT46mnnlrrfZ5++umYOHFip2377bdfXH/99dHS0hI1NTVr3KepqSmamprarzc2NvZA9QBATyp1v66uyr43gvkfi3yvammLxStaYvHK5li8oiWWN7Vu9LQuFgftvoZcVQyor40BDTUxsKE26muNKgcoF8Xo2XU1VbHVoIbYalBD5POFWLyyJRataI7FK5pjycoWcyQnQFVVJgbUv9fHBzbURr/6amuQEDVV/g9UiuqqTNF/5lMVnM+fPz/a2tpi2LBhnbYPGzYs5s6du9b7zJ07d637t7a2xvz582PEiBFr3Gfq1Klx/vnn91zhAECPK7d+XVdTFcP7V8Xw/u8F6a1t+VjW1Bp53yztso+MHtjt+2Qz7402rKkyCg2gXBW7Z2ezmdi8T21s3qc2IiLy+UI0rmqJRSveC9OXrGiJth5o0NYk2TTVVZn2kHxAn5romxOUs6Zqf+NVjFL8PZ+q4Hy19/8iLRQK6/3lurb917Z9tSlTpsRZZ53Vfr2xsTFGjRq1seUCAL2g3Pt1dVU2BjTUbtR9/3zBfj1cTTI01KbyT1eAilfqnp3NZmJAQ20MaKiNsdGnPUhfPd3akpUtG/UtsVKvB5K0NUmy2YgBDbUxqE9tDOxTKyinS4w4rxyC8000ePDgqKqqWmN0+bx589YYVb7a8OHD17p/dXV1DBo0aK33yeVykcvleqZoAKBXpLlfC5ABSJNy69kdg/Sth0S0tOVj0fLmmL+sORYsb4qmFvO69JSG2qoY3DcXg/q893oXe+E/kq+UC+Bad6i4SnGSJFWfumpra2O33XaLhx9+OA499ND27Q8//HD8+7//+1rvs9dee8V9993Xaduvf/3r2H333dc6vzkAAABQOWqqsjG0X10M7VcXhUIhlja1xrtLm2LBsuZoXNmyzvtVaqgWse5gLZOJGNBQE4M3y8XgzXLRJ5eqWIoSqK0q3bcqrDtUXKWYeip1v6HOOuus+MIXvhC777577LXXXnHdddfFzJkz46STToqI974CNmvWrLjpppsiIuKkk06KH/7wh3HWWWfFiSeeGE8//XRcf/31cfvtt5fyaQAAAABlJpPJRL+6muhXVxPbDHkvnF6wrDlWriUk/8SlT5Sgwn8qVagWEfG/X9tnjW2b5apj8z611h2hR40YUBeDNtu46Q/ZOHtts/YZOnpbtRHnm+6oo46KBQsWxAUXXBBz5syJcePGxYMPPhijR4+OiIg5c+bEzJkz2/cfO3ZsPPjgg3HmmWfGVVddFSNHjowrrrgiDj/88FI9BQAAoIKtaG4tdQklYRoqkihXXRUjB9SXuoyys82QzUpdAhWipipbspMx1h1Kv0yhsDFLXNBRY2Nj9O/fP5YsWRL9+vUrdTkApIw+0zO8jkBvKxQKsaypNZauao1N+ZD1r997tMdqSpLf/ee+G33f2qpsDGioKflIVr2mZ6TldazUk2ARlRWsAcnT1T7jNxlAhajUP9z90Q5Ab1o9TcPC5c2xYHlztLRaNHBjvTa7cZPun8lE9Kuvic371MbgPrnoV18dmYyFBikdf4cCJFuXfosfdthhXX7AX/ziFxtdDEAlyOcLkS/Bl322//ZDRT9mOXjjwgNKctyqbMaHdYCUWtHcGvMam2Le0qb1Lgy4sa46etcef8wNaWpta58P+bIjdirJAlybqlCIWLKiJZasaIkZ7y6PmupsDN6sNob2rYtBfWojm9WXAYCu61Jw3r9//96uA6BizFy4It6Yt6zUZVSMx6e/W5LjfnTM5tG/oaYkxwag57W25WNu46qYs2RVLFnR82F5R7ma0obWueqqktfQE1pa8zFn8aqYs3hVVFdlYnj/uthiQH30rdOfAYAN61JwPm3atN6uA6BiNLeV5ivcRq8VV1NbW0T4YA6QdM2t+Zi5cHm8tWhltLVZHiqpWtsK8fbClfH2wpUxsE9NjBnUJwZtlit1WQBAGdvoCbfefffdmD59emQymfjQhz4UQ4YM6cm6AFKruURzn5Z65FhaRq91VaneZwB6zlsLV8Qb7y4remDe1NJW1ONFvHeye23/LqZi/Z2waHlLLFq+ODbfrDa2H9Ev6iro7xMAoOu6HZwvX748TjvttLjpppsin38vFKiqqopjjjkmrrzyymhoaOjxIgHSpFQjzimuFqMSARLt1dlLYs7iVSU59lduf7Ekx11t9bfFiu1/jtm9qMdbuKw5fj9jYXxkqwGmbwEA1pDt7h3OOuuseOKJJ+K+++6LxYsXx+LFi+Oee+6JJ554Ir72ta/1Ro0AqdJiJHJFaHGCBCCx5i1dVbLQnOJqac3HX+YuLXUZAEAZ6vaI87vuuit+/vOfx/jx49u3ffrTn476+vo48sgj4+qrr+7J+gBSx0jkymCqFoDkKpS4VVuXpLjyeX+bAQBr6nZwvmLFihg2bNga24cOHRorVqzokaIA0sxI5MrgfQZIrqF9czG8f13MXVKaUeelXhOkktYlqa3OxnYj+5W6DACgDHV7qpa99torzj333Fi16p9/RK5cuTLOP//82GuvvXq0OIC0yecL0WZUU0XwzQKA5MpkMjFui/6x05b9o0+u22ONSIBsNmLLzetjz60HRT/zmwMAa9HtvwL/3//7f7H//vvHlltuGTvvvHNkMpl46aWXoq6uLh566KHeqBEgNVryRiFXilYjzgESb2i/uhjary4WLm+O2YtXxrvLmqKtCCdGm1raev0YaxyztW2t/y6mYoxy71tXHSP618fw/nVRW93tcWQAQAXpdnA+bty4eP311+OWW26Jv/zlL1EoFGLSpEnx+c9/Purr63ujRoDUaDUKuWK0+GYBQGps3qc2Nu9TG/l8IRataI75y5pjwfKmWNHUOwHzV25/sVcet6tWz3VebP9zzO49/phVVZkY2FAbg/rUxpC+uairkCloAIBNt1HfO6yvr48TTzyxp2sBSD3BeeUw4hwgfbLZTAzaLBeDNstFRN9Y1dIWS1a2xIrm0ozQTptthm7WY49VlclE/4aa6FdXHZlMpsceFwCoHBsVnM+aNSt+97vfxbx58yL/vmkHTj/99B4pDCCNWk3VUjEKhYi2fCGqsj6sA6RVXU1Vr4xg/vMF+/X4YyZBQ6355AGA8tHtv0ymTZsWJ510UtTW1sagQYM6nb3PZDKCc4D1aDV9R0VpactHVdZXwgHoHgEyAEDpdfsvsm9/+9vx7W9/O6ZMmRLZrMVUALpDcF5Z8gXvNwAAACRRt5PvFStWxKRJk4TmABuhzRznFcWJEgAAAEimbqffxx9/fPzsZz/rjVoAUq/NCOSK4kQJAAAAJFO3p2qZOnVqHHjggfGrX/0qdtxxx6ipqel0+2WXXdZjxQGkTZsRyBXFiRIAAABIpm4H59/97nfjoYceig9/+MMREWssDgrAupnzurLknSgBAACAROp2cH7ZZZfFDTfcEJMnT+6FcgDSzYjzyuLtBgAAgGTq9hznuVwu/vVf/7U3agFIPcF5ZTFVCwAAACRTt4PzM844I6688sreqAUg9eSolcVULQAAAJBM3Z6q5dlnn41HH3007r///thhhx3WWBz0F7/4RY8VB5A25jivLN5uAAAASKZuB+cDBgyIww47rDdqAUg9wXll8X4DAABAMnU5OF+2bFlsttlmMW3atN6sByDVxKiVxfsNAAAAydTlOc4HDx4cBxxwQFx99dUxe/bs3qwJILUKRiBXFCPOAQAAIJm6HJxPnz49Pv3pT8ddd90VY8eOjY9+9KPx3//93/Hyyy/3Zn0AqSJHrSzebwAAAEimLgfno0ePjtNOOy0eeeSRmDdvXpx11lnx6quvxsc//vEYO3ZsnHHGGfHoo49GW1tbb9YLAAkiOQcAAIAk6nJw3lH//v3j6KOPjjvuuCPmz58f1157beTz+TjuuONiyJAhceutt/Z0nQCpIEatLEacAwAAQDJ1eXHQdT5AdXVMnDgxJk6cGFdeeWW8+OKL0dra2hO1AaSOIBUAAACg/HU5OO/KXObV1dUxevTo2HzzzTepKAAAAAAAKJUuB+e77LJLZDKZKGxguGQmk4lddtklfvKTn8S4ceM2uUAAAAAAACimLgfnM2bM2OA++Xw+3nnnnbj44ovj5JNPjieffHKTigNIm0ym1BUAAAAAsCFdDs5Hjx4dERHNzc1RW1u71n3mz58fY8eOjYsuuih23nnnnqkQABLKiRIAAABIpmx373DkkUdGPp9fY/s777wT48ePj4iIsWPHxlNPPbXJxQGkjRwVAAAAoPx1OzifM2dOHH/88WtsGz9+fGy77bYREVFVVWXEOcBaGIFcabzhAAAAkETdDs4ffPDBePbZZ+PMM8+MiIhZs2bF+PHjY8cdd4yf/vSnPV4gQLoIUiuJEyUAAACQTF2e43y1QYMGxUMPPRQf+9jHIiLigQceiI985CNx6623Rjbb7RweoKIIUitL1hsOAAAAidTt4DwiYsstt4yHH344Pvaxj8WnPvWpuPnmmyMjHADYIEFqZcl6uwEAACCRuhScDxw4cK3B+IoVK+K+++6LQYMGtW9buHBhz1UHkDKC1MqSMTUPAAAAJFKXgvPLL7+8l8sAqAxGnFcWM5gBAABAMnUpOD/22GN7uw6AiiA3ryxOlAAAAEAyGQsHUERV5mqpKIJzAAAASCbBOUARCVIri6laAAAAIJl8pAcoIsF5ZanyfgMAAEAiCc4BishULZXF+w0AAADJ1OXgfOTIkXHyySfHL3/5y2hubu7NmgBSq1qQWlEE5wAAAJBMXQ7Ob7vttmhoaIjTTz89Bg8eHEcccUTcfPPNsXDhwt6sDyBVsoLUiiI4BwAAgGTqcnA+fvz4uPTSS+P111+Pp59+Oj7ykY/EVVddFSNGjIjx48fHD37wg/i///u/3qwVIPGMOK8sgnMAAABIpo2a43yHHXaIKVOmxDPPPBMzZ86Mz3/+8/Hoo4/GjjvuGOPGjYsHHnigp+vcoDfffDOOP/74GDt2bNTX18c222wT55577ganlZk8eXJkMplOlz333LNIVQOVRpBaWaqzlhIBAACAJKre1AcYNmxYnHjiiXHiiSfGihUr4qGHHopcLtcTtXXLX/7yl8jn83HttdfGBz7wgXjllVfixBNPjOXLl8cll1yy3vvuv//+MW3atPbrtbW1vV0uUKGMOK8cmYwTJQAAAJBUmxycd9TQ0BCHHnpoTz5kl+2///6x//77t1/feuutY/r06XH11VdvMDjP5XIxfPjw3i4RIKqrjECuFEJzAAAASK4eDc7LzZIlS2LzzTff4H6PP/54DB06NAYMGBD77LNPXHjhhTF06NB17t/U1BRNTU3t1xsbG3ukXiD9jDivHKZpKT39GgCSQc8GoByl9lP9//3f/8WVV14ZJ5100nr3O+CAA+LWW2+NRx99NC699NJ47rnnYt999+3UtN9v6tSp0b9///bLqFGjerp8IKUE55Wjusp7XWr6NQAkg54NQDkq++D8vPPOW2Pxzvdfnn/++U73mT17duy///5xxBFHxAknnLDexz/qqKPiM5/5TIwbNy4OOuig+OUvfxl//etf17vA6ZQpU2LJkiXtl7feeqtHniuQftVV2cjIUytCjeC85PRrAEgGPRuActStqVpaW1ujrq4uXnrppRg3blxv1dTJqaeeGpMmTVrvPmPGjGn/9+zZs2PChAmx1157xXXXXdft440YMSJGjx4dr7/++jr3yeVyJVkAFUiHqmwmWtsKpS6DXmaqltLTrwEgGfRsAMpRt4Lz6urqGD16dLS1tfVWPWsYPHhwDB48uEv7zpo1KyZMmBC77bZbTJs2LbIbEVosWLAg3nrrrRgxYkS37wvQFbVV2Wgt4u9RSsNULQAAAJBc3U6Wv/nNb8aUKVNi4cKFvVHPRps9e3aMHz8+Ro0aFZdcckm8++67MXfu3Jg7d26n/bbddtu4++67IyJi2bJlcfbZZ8fTTz8db775Zjz++ONx0EEHxeDBg+PQQw8txdMAKkB1lZHIlaDW+wwAAACJ1a0R5xERV1xxRbzxxhsxcuTIGD16dPTp06fT7X/4wx96rLju+PWvfx1vvPFGvPHGG7Hlllt2uq1Q+OeUCNOnT48lS5ZERERVVVX86U9/iptuuikWL14cI0aMiAkTJsSdd94Zffv2LWr9QOUwErkyOEECAAAAydXt4PyQQw7phTI23eTJk2Py5Mkb3K9jiF5fXx8PPfRQL1YFsCYjkSuDxUEBAAAgubodnJ977rm9UQdAxaitFpxXAidIAAAAILm6HZyv9sILL8Rrr70WmUwmtt9++9h11117si6A1KoRqFYE7zMAAAAkV7eD83nz5sWkSZPi8ccfjwEDBkShUIglS5bEhAkT4o477oghQ4b0Rp0AqWEKj8pQ45sFAAAAkFjd/lR/2mmnRWNjY7z66quxcOHCWLRoUbzyyivR2NgYp59+em/UCJAqpvCoDN5nAAAASK5ujzj/1a9+FY888khst9127du23377uOqqq2LixIk9WhxAGpnjPP0yGd8sAAAAgCTrdnqTz+ejpqZmje01NTWRz+d7pCiANBOcp19NVTYyGcE5AAAAJFW305t99903zjjjjJg9e3b7tlmzZsWZZ54Zn/jEJ3q0OIA0smhk+jk5AgAAAMnW7U/2P/zhD2Pp0qUxZsyY2GabbeIDH/hAjB07NpYuXRpXXnllb9QIkCo1VdnIylVTTXAOAAAAydbtOc5HjRoVf/jDH+Lhhx+Ov/zlL1EoFGL77bePT37yk71RH0Aq1VZVxap8W6nLoJdYGBQAAACSrVvBeWtra9TV1cVLL70Un/rUp+JTn/pUb9UFkGq5mmysahGcp1XOiHMAAABItG59sq+uro7Ro0dHW5uwB2BTGJGcbrnqqlKXAAAAAGyCbic33/zmN2PKlCmxcOHC3qgHoCLkagTnaWaOcwAAAEi2bs9xfsUVV8Qbb7wRI0eOjNGjR0efPn063f6HP/yhx4oDSCsjktPNVC0AAACQbN0Ozg855JBeKAOgshiRnG6+UQAAAADJ1u3FQSMivvjFL8aoUaN6pSCASlAnOE813ygAAACAZOv24qCXXHKJxUEBNlGuRrCaVtVVmajKZkpdBgAAALAJuj3k8ROf+EQ8/vjjvVAKQOUwB3Z6GW0OAAAAydftOc4POOCAmDJlSrzyyiux2267rbE46MEHH9xjxQGkVU1VNqqymWjLF0pdCj3M/OYAAACQfN0Ozk8++eSIiLjsssvWuC2TyZjGBaCLctXZWNHsd2ba1BlxDgAAAInX7eA8n8/3Rh0AFSdXUyU4T6E6I84BAAAg8Xy6BygRAWs61Vn4FQAAABKvy6nNpz/96ViyZEn79QsvvDAWL17cfn3BggWx/fbb92hxAGkmYE0n7ysAAAAkX5eD84ceeiiamprar1900UWxcOHC9uutra0xffr0nq0OIMUErOnkmwQAAACQfF3+dF8oFNZ7HYDuqasWsKZRzuKgAAAAkHhSG4ASMeI8fWqrs1GVzZS6DAAAAGATdTk4z2Qykclk1tgGwMYRnKeP9xQAAADSobqrOxYKhZg8eXLkcrmIiFi1alWcdNJJ0adPn4iITvOfA7BhVdlM1FZno7k1X+pS6CHmNwcAAIB06HJwfuyxx3a6/h//8R9r7HPMMcdsekUAFaSupkpwniL1RpwDAABAKnQ5OJ82bVpv1gFQkeprqqJxZUupy6CHmKoFAAAA0qHLwTkAPa++tnhTezS1tBXtWO3HbG1b67+LKVfEMFtwDgAAAOkgOAcooVx18YLWr9z+YtGOtTZn/ezlkhz3f47ZvWjHqq8VnAMAAEAaWMUMoIQErelSV62tAgAAQBoYcQ5QQsVcTPKqo3ct2rFWa2ptax9pftkROxV1hH2x1VRno7pKcA4AAABpIDgHKKFiBufFnOt7rcevrip5Db2pmO8lAAAA0LsMjQMooWw2E7kav4rToMG0OwAAAJAa0hqAEjNSOR3qvI8AAACQGoJzgBITuKaDhV4BAAAgPQTnACVmio908M0BAAAASA/BOUCJGamcDk6AAAAAQHoIzgFKzEjl5MtmI3LVWioAAACkhU/5ACVmxHny1VVXRSaTKXUZAAAAQA8RnAOUWK66KqqyQtckc/IDAAAA0kVwDlAGBK/J5v0DAACAdBGcA5QBC0smW0NNdalLAAAAAHqQ4BygDFggNNmMOAcAAIB0EZwDlAHBa7J5/wAAACBdBOcAZaCh1lQfSdbgGwMAAACQKoJzgDJgjvPkqqupimw2U+oyAAAAgB4kOAcoA7nqbGT9Rk4k07QAAABA+ohpAMpAJpOJOtN9JJJvCwAAAED6CM4BykQf85wnkuAcAAAA0idVwfmYMWMik8l0uvznf/7neu9TKBTivPPOi5EjR0Z9fX2MHz8+Xn311SJVDPBPpvxIJu8bAAAApE+qgvOIiAsuuCDmzJnTfvnmN7+53v2///3vx2WXXRY//OEP47nnnovhw4fHpz71qVi6dGmRKgZ4j5HLydTgmwIAAACQOqkLzvv27RvDhw9vv2y22Wbr3LdQKMTll18e3/jGN+Kwww6LcePGxU9+8pNYsWJF3HbbbUWsGiCi3hznidTgfQMAAIDUSV1wftFFF8WgQYNil112iQsvvDCam5vXue+MGTNi7ty5MXHixPZtuVwu9tlnn3jqqafWeb+mpqZobGzsdAHYVH1yRi4nTV1NVWSzmVKXwTro1wCQDHo2AOUoVcH5GWecEXfccUc89thjceqpp8bll18ep5xyyjr3nzt3bkREDBs2rNP2YcOGtd+2NlOnTo3+/fu3X0aNGtUzTwCoaLnqbGRT9Vs5/cxvXt70awBIBj0bgHJU9hHNeeedt8aCn++/PP/88xERceaZZ8Y+++wTO+20U5xwwglxzTXXxPXXXx8LFixY7zEymc6jBQuFwhrbOpoyZUosWbKk/fLWW29t+hMFKl4mk4n6GqPOk8S89OVNvwaAZNCzAShHZZ/QnHrqqTFp0qT17jNmzJi1bt9zzz0jIuKNN96IQYMGrXH78OHDI+K9kecjRoxo3z5v3rw1RqF3lMvlIpfLbah0gG5rqK2K5U2tpS6DLhKclzf9GgCSQc8GoByVfXA+ePDgGDx48Ebd98UXX4yI6BSKdzR27NgYPnx4PPzww7HrrrtGRERzc3M88cQTcdFFF21cwQCbQBCbLKZqAQAAgHQq+6lauurpp5+OH/zgB/HSSy/FjBkz4qc//Wl8+ctfjoMPPji22mqr9v223XbbuPvuuyPivWkRvvrVr8Z3v/vduPvuu+OVV16JyZMnR0NDQ3zuc58r1VMBKpggNln61Jb9+WcAAABgI6TmE38ul4s777wzzj///GhqaorRo0fHiSeeGOecc06n/aZPnx5Llixpv37OOefEypUr45RTTolFixbFHnvsEb/+9a+jb9++xX4KANEgiE2U+honOgAAACCNUpPQfOQjH4lnnnlmg/sVCoVO1zOZTJx33nlx3nnn9VJlAF1nqpbkqK+timx23QtJAwAAAMmVmqlaANKgrqYqqoSxiWBaHQAAAEgvwTlAmRHIJoNvBwAAAEB6Cc4ByoxANhksDAoAAADpJTgHKDMWCE0G3wwAAACA9BKcA5QZI86TwfsEAAAA6SU4BygzAtnyl8lE1Nd4nwAAACCtBOcAZcYUIOWvvqYqMplMqcsAAAAAeongHKDM5KqroqpKKFvOGnLmoQcAAIA0E5wDlKEG04CUNdPpAAAAQLoJzgHKUB8jmsua+c0BAAAg3QTnAGXIPOflzYhzAAAASDfBOUAZEsyWN98IAAAAgHQTnAOUoYYawWy5ymYjctXaJwAAAKSZT/4AZaghZ8R5uaqvqY5MJlPqMgAAAIBeJDgHKEM1VdmorhLOliPT6AAAAED6Cc4BypR5tMuT4BwAAADST3AOUKbqawS05ajBCQ0AAABIPcE5QJkysrk8NTihAQAAAKknOAcoUw21RjaXo3onNAAAACD1BOcAZUpAW36qspmoM+IcAAAAUk9wDlCmTNVSfpzMAAAAgMogOAcoUzVV2aip9mu6nDiZAQAAAJVBIgNQxgS15cX7AQAAAJVBcA5QxurNp11W6i3YCgAAABVBcA5QxvrkBLXlpMGJDAAAAKgIgnOAMmbEeXmxOCgAAABUBsE5QBkT1JaPqmwm6pzIAAAAgIogOAcoYxajLB9OYgAAAEDlEJwDlLGaqmzUVPtVXQ5MmwMAAACVQxoDUOaMOi8P3gcAAACoHIJzgDJnpHN5MFULAAAAVA7BOUCZE9iWh4ba6lKXAAAAABSJ4BygzJkipDx4HwAAAKByCM4BylxDjZHOpZbNRuQs0goAAAAVQwoAUOZM1VJ6dTVVkclkSl0GAAAAUCSCc4AyV1udjaoqoW0pmd8cAAAAKovgHCAB6muMOi8lrz8AAABUFsE5QAJYmLK0vP4AAABQWQTnAAkguC2tOiPOAQAAoKIIzgESQHBbWk5cAAAAQGURnAMkgMUpS8sc5wAAAFBZBOcACSC4LZ1cTTay2UypywAAAACKSHAOkAB1NdnIyG5LwjQtAAAAUHkE5wAJkMlkjDovEfPLAwAAQOURnAMkRJ2RzyVhfnkAAACoPIJzgIQw4rw0vO4AAABQeQTnAAlhru3SEJwDAABA5RGcAySEALc06p2wAAAAgIojOAdICAFu8VVVZaK2WqsEAACASiMNAEgII86Lz2sOAAAAlSk1wfnjjz8emUxmrZfnnntunfebPHnyGvvvueeeRawcoGuqq7JRXZUpdRkVRXAOAAAAlam61AX0lL333jvmzJnTadu3vvWteOSRR2L33Xdf733333//mDZtWvv12traXqkRYFM11FZH48qWUpdRMUyPAwAAAJUpNcF5bW1tDB8+vP16S0tL3HvvvXHqqadGJrP+EZq5XK7TfQHKVX1NleC8iIw4BwAAgMqUmuD8/e69996YP39+TJ48eYP7Pv744zF06NAYMGBA7LPPPnHhhRfG0KFD17l/U1NTNDU1tV9vbGzsiZIBNqi+NjUzbCVCneA80fRrAEgGPRuAcpTaBOb666+P/fbbL0aNGrXe/Q444IC49dZb49FHH41LL700nnvuudh33307Ne33mzp1avTv37/9sqFjAPQUQW5xmaol2fRrAEgGPRuAclT2wfl55523zkU/V1+ef/75Tvd5++2346GHHorjjz9+g49/1FFHxWc+85kYN25cHHTQQfHLX/4y/vrXv8YDDzywzvtMmTIllixZ0n556623Nvl5AnSFqUOKy+udbPo1ACSDng1AOSr7qVpOPfXUmDRp0nr3GTNmTKfr06ZNi0GDBsXBBx/c7eONGDEiRo8eHa+//vo698nlcpHL5br92ACbygjo4qmtzkZVdv1rZFDe9GsASAY9G4ByVPbB+eDBg2Pw4MFd3r9QKMS0adPimGOOiZqamm4fb8GCBfHWW2/FiBEjun1fgN5WVy04LxYnKQAAAKBylf1ULd316KOPxowZM9Y5Tcu2224bd999d0RELFu2LM4+++x4+umn480334zHH388DjrooBg8eHAceuihxSwboEuy2UzkalL3q7ssmaYFAAAAKlfZjzjvruuvvz723nvv2G677dZ6+/Tp02PJkiUREVFVVRV/+tOf4qabborFixfHiBEjYsKECXHnnXdG3759i1k2QJfV11RFU0u+1GWkXp0TFAAAAFCxUhec33bbbeu9vVAotP+7vr4+Hnrood4uCaBH1dVURURLqctIvTojzgEAAKBiGU4HkDAC3eIwVQsAAABULsE5QMJYtLI4vM4AAABQuQTnAAljJHRx1FV7nQEAAKBSCc4BEkZw3vtyNdnIZjOlLgMAAAAoEcE5QMLkqrORken2KicnAAAAoLIJzgESJpvNRM40Ir3KAqwAAABQ2QTnAAlUV+PXd28SnAMAAEBlk7wAJJBgt3fV13p9AQAAoJIJzgESSHDeu+qqtUcAAACoZJIBgAQyIrp3eX0BAACgsgnOARLIiOjeVWfxVQAAAKhokheABDIiuvfUVmcjm82UugwAAACghATnAAlkRHTvcVICAAAAEJwDJFA2m4la07X0CiclAAAAAKkLQEIZGd076mq0RgAAAKh00gGAhDIyunfU1XhdAQAAoNIJzgESysjo3iE4BwAAAKQuAAkl4O0dTkgAAAAA0gGAhBKc9456rysAAABUPME5QEIZGd3zqqsyUV3ldQUAAIBKJx0ASCgjo3ueUfwAAABAhOAcILGqq7JRXZUpdRmpIjgHAAAAIgTnAIkm6O1ZRvEDAAAAEYJzgEQTnPcs88YDAAAAEYJzgEQT9PYsJyIAAACACME5QKKZWqRn1VV7PQEAAADBOUCiGSHds3JG8AMAAAAhOAdINCOke042G5Gr1hYBAAAAwTlAohkh3XPqqqsik8mUugwAAACgDEhcABIsV52NrN/kPSJn2hsAAADgH8QtAAmWyWQiZ7qWHlFn9D4AAADwD1ICgIQT+PYMC60CAAAAq0lbABLOiPOeITgHAAAAVhOcAyRcfa3AtyfUVWuJAAAAwHukBAAJZ6R0z/A6AgAAAKsJzgESzkjpniE4BwAAAFaTtgAknMB309VUZ6Mqmyl1GQAAAECZEJwDJJzgfNMZtQ8AAAB0JCkASLiqbCZqBL+bxAKrAAAAQEeSFoAUMGJ60xi1DwAAAHQkaQFIAcHvpqmr9voBAAAA/yQ4B0gBU41smroa7RAAAAD4J0kBQAoYMb1pckbsAwAAAB0IzgFSwIjpTVMvOAcAAAA6kLQApIAR0xsvm42otbgqAAAA0IGkACAFjDjfeKa5AQAAAN5P0gKQArnqqsj6jb5R6iysCgAAALyPmAUgJYyc3jheNwAAAOD9BOcAKWHk9MYxzQ0AAADwftICgJQwcnrj1FlYFQAAAHgfwTlAShg5vXEE5wAAAMD7JSZlufDCC2PvvfeOhoaGGDBgwFr3mTlzZhx00EHRp0+fGDx4cJx++unR3Ny83sdtamqK0047LQYPHhx9+vSJgw8+ON5+++1eeAYAvUsAvHHqvW4AAADA+yQmOG9ubo4jjjgiTj755LXe3tbWFp/5zGdi+fLl8dvf/jbuuOOOuOuuu+JrX/vaeh/3q1/9atx9991xxx13xG9/+9tYtmxZHHjggdHW1tYbTwOg1wjON06uOjGtEAAAACiS6lIX0FXnn39+RETceOONa73917/+dfz5z3+Ot956K0aOHBkREZdeemlMnjw5LrzwwujXr98a91myZElcf/31cfPNN8cnP/nJiIi45ZZbYtSoUfHII4/Efvvt1ztPBqAXGDndfbmabGSzmVKXAQAAAJSZ1Ayze/rpp2PcuHHtoXlExH777RdNTU3xwgsvrPU+L7zwQrS0tMTEiRPbt40cOTLGjRsXTz311DqP1dTUFI2NjZ0uAKVm5HT3OdmQbvo1ACSDng1AOUpNyjJ37twYNmxYp20DBw6M2tramDt37jrvU1tbGwMHDuy0fdiwYeu8T0TE1KlTo3///u2XUaNGbfoTANhE2WwmchYI7RbT26Sbfg0AyaBnA1COSpqwnHfeeZHJZNZ7ef7557v8eJnMml+3LxQKa92+Phu6z5QpU2LJkiXtl7feeqtbjw/QW4yg7p46JxpSTb8GgGTQswEoRyWd4/zUU0+NSZMmrXefMWPGdOmxhg8fHr///e87bVu0aFG0tLSsMRK9432am5tj0aJFnUadz5s3L/bee+91HiuXy0Uul+tSXQDF9N4I6pZSl5EYuWonGtJMvwaAZNCzAShHJQ3OBw8eHIMHD+6Rx9prr73iwgsvjDlz5sSIESMi4r0FQ3O5XOy2225rvc9uu+0WNTU18fDDD8eRRx4ZERFz5syJV155Jb7//e/3SF0AxWQEdffU1wrOAQAAgDUlJmGZOXNmvPTSSzFz5sxoa2uLl156KV566aVYtmxZRERMnDgxtt9++/jCF74QL774Yvzv//5vnH322XHiiSdGv379IiJi1qxZse2228azzz4bERH9+/eP448/Pr72ta/F//7v/8aLL74Y//Ef/xE77rhjfPKTnyzZcwXYWObs7h6vFwAAALA2JR1x3h3f/va34yc/+Un79V133TUiIh577LEYP358VFVVxQMPPBCnnHJK/Ou//mvU19fH5z73ubjkkkva79PS0hLTp0+PFStWtG/7wQ9+ENXV1XHkkUfGypUr4xOf+ETceOONUVUlTAGSRxDcPeaEBwAAANYmMcH5jTfeGDfeeON699lqq63i/vvvX+ftY8aMiUKh0GlbXV1dXHnllXHllVf2RJkAJSU477qa6mxUZbu3eDQAAABQGRIzVQsAG2YEddfVVWuBAAAAwNpJDQBSpCqbiRqBcJdYGBQAAABYF+kKQMoYdd41prUBAAAA1kVwDpAydTV+tXeFEwwAAADAukhXAFJGINw1OScYAAAAgHWQGgCkjClIusYJBgAAAGBdqktdAAA9a13BeVNLW5EriWhqbVvrv4spt47XwwkGAAAAYF0E5wApU1+79kD4K7e/WORKOjvrZy+X5Lj/c8zua2yrrspETZUvXQEAAABrJzgHSJk+tVWx1zaDSl1G2Vjba5HNZEpQCQAAAJAUgnOAlMlkMtEnt+av9z9fsF8Jqim9hlqtDgAAAOgeaQJAhRAgAwAAAHSNCV4BAAAAAKADwTkAAAAAAHQgOAcAAAAAgA4E5wAAAAAA0IHgHAAAAAAAOhCcAwAAAABAB4JzAAAAAADoQHAOAAAAAAAdCM4BAAAAAKADwTkAAAAAAHQgOAcAAAAAgA4E5wAAAAAA0IHgHAAAAAAAOhCcAwAAAABAB4JzAAAAAADoQHAOAAAAAAAdCM4BAAAAAKCD6lIXkAaFQiEiIhobG0tcCQBptLq/rO43bBz9GoDepmf3DD0bgN7U1X4tOO8BS5cujYiIUaNGlbgSANJs6dKl0b9//1KXkVj6NQDFomdvGj0bgGLYUL/OFJwK32T5fD5mz54dffv2jUwmU+pyEqmxsTFGjRoVb731VvTr16/U5QA9yM/3pisUCrF06dIYOXJkZLNmWdtY+vWm8/MM6eZnfNPp2T1Dz940fpYh3fyMb7qu9msjzntANpuNLbfcstRlpEK/fv380ENK+fneNEatbTr9uuf4eYZ08zO+afTsTadn9ww/y5BufsY3TVf6tVPgAAAAAADQgeAcAAAAAAA6EJxTFnK5XJx77rmRy+VKXQrQw/x8Q3r4eYZ08zMO6eBnGdLNz3jxWBwUAAAAAAA6MOIcAAAAAAA6EJwDAAAAAEAHgnMAAAAAAOhAcA4AAAAAAB1Ul7oAOO644za4T6FQiBtvvLH3iwF6XFd+xiMipk2b1suVAJtCv4Z0068hPfRsSC/9urgyhUKhUOoiqGyHHXbYOm8rFArx4osvxsyZMyOfzxexKqCnVFVVxT777BMDBgxY736/+MUvilMQsFH0a0g3/RrSQ8+G9NKvi0twTllasGBB3HLLLTFt2rSYOXNmTJo0KX70ox+VuixgI1RVVcUf/vCH2HnnnUtdCtDD9GtID/0a0k3PhnTQr4vLHOeUjXw+Hw8++GAcccQRseWWW8b9998fX//612POnDkaOgCUCf0aAJJBzwbYNOY4pywsW7Ystt1226itrY3JkyfHpZdeGltttVWpywIAOtCvASAZ9GyATWfEOWUjm81GJpOJQqEQZhACgPKkXwNAMujZAJtGcE5Z2GyzzeLvf/97/OhHP4pXX301tt9++/jUpz4Vt99+ezQ1NZW6PGATHHvssTFw4MBSlwH0AP0a0ku/hnTRsyGd9OvisjgoZWnhwoVxyy23xA033GDhEki4rbfeeoMjXAqFQrz55pvFKQjoMfo1pId+DemmZ0M66NfFJTin5A499NB13pbJZOLFF1+Mv//975HP54tYFdBTrrjiinXe9s4778Stt94aM2fO9DMOZU6/hnTTryE99GxIL/26uCwOSsltvvnm6z1bNmHChCJWA/S0008/vdP1tra2eOCBB2LatGnx0EMPxT777BPf//73S1Qd0FX6NaSbfg3poWdDeunXxWXEOQBF8Ze//CVuuOGGuOWWW6JPnz5x3HHHxbHHHhtbbLFFqUsDAP5BvwaA8qdfF4cR55Tccccdt8F9CoVC3Hjjjb1fDNArFi1aFDvssEN8/OMfjzvuuCM+/vGPl7okoJv0a0g//RrSQc+GdNOviydb6gJgyZIl7ZfZs2fHLbfc0mnbvHnz4qabbip1mcAmaGhoiEmTJsWzzz4b5557btx0002xcuXKUpcFdIN+DemnX0M66NmQbvp18ZiqhbIyY8aM2GmnnWLp0qXt2959990YPnx4tLW1lbAyoCc0NjbGrbfeGtdff3288cYbcdRRR8UXv/jF2GOPPUpdGtAN+jWkm34N6aFnQ3rp173PiHPKSnV1dbS2tnba1tzcHFVVVSWqCOhJ/fr1i5NPPjmef/75+M1vfhP19fXxmc98JnbYYYdSlwZ0g34N6aZfQ3ro2ZBe+nXvM8c5ZWXkyJHR2toazz77bPzLv/xLREQ8+eSTMWrUqBJXBmyszTffPNb15aZMJhOtra3x2muvFbkqYFPo15A++jWkk54N6aJfF5fgnLJSVVUVhxxySHzmM5+Jz33uc7Fy5cq49dZb45RTTil1acBGuvzyy0tdAtDD9GtIH/0a0knPhnTRr4vLHOeUnfnz58fpp58ejzzySNTU1MQhhxwSF198cTQ0NJS6NADgH/RrAEgGPRtg4wjOAQAAAACgA1O1UHLnn39+fPazn13n4gXTpk2L3XbbLXbaaaciVwb0hK233nqdc7CtVigU4s033yxOQcBG0a8h3fRrSA89G9JLvy4uI84puaqqqhg0aFA89thja23sX/3qV2PevHlx2223laA6YFNdccUV7f+eP39+XHzxxXHRRRe1b1u2bFl885vfjHw+X4rygC7SryHd9GtIDz0b0ku/Li7BOSVXVVUVn//85+PXv/51PProo7H99tt3uv3JJ5+Mz3/+8zFz5swSVQj0lL/97W+x8847x9KlS9u3zZs3L4YPH66xQ5nTr6Fy6NeQbHo2VAb9uvdlS10ARERcfPHFcdxxx8W+++4br776aqfbRo0aFfPnzy9RZUBPamhoiKampmhra2vftmzZsqivry9hVUBX6ddQGfRrSD49G9JPv+59gnPKxtSpU+NLX/pS7LvvvvH000+3b3/qqadi7NixJawM6CnDhw+PXC4XN998c/u2adOmxdZbb13CqoDu0K8h/fRrSAc9G9JNv+59FgelrFxwwQVRX18fEyZMiCOOOCL69+8f06ZNi6lTp5a6NKCHnHrqqXH88cfHD37wg1i1alW8/vrrcc0115S6LKAb9GtIP/0a0kHPhnTTr3uXOc4pua233jqeffbZGDx4cPu2Rx99NK655ppYuHBhHHbYYXHKKaeUsEKgp1133XXxyCOPRE1NTRxyyCFxxBFHlLokYAP0a6g8+jUkk54NlUW/7j2CcwAAAAAA6MAc5wD0qgkTJsTtt9++ztvPPPPMuPvuu4tYEQDwfvo1AJQ//bq4BOcA9Krf/OY3ceyxx66zuW+xxRZx9dVXF7kqAKAj/RoAyp9+XVyCcwB63RVXXBEnnHBC3HnnnWvcNmHChHjppZeKXxQA0Il+DQDlT78unupSFwBA+h166KGxxRZbxNFHHx35fD6OPvro9tsGDBgQTU1NJawOAIjQrwEgCfTr4jHiHICiOOigg+LnP/95fOlLX4qLL744Vq9NfcMNN8QOO+xQ4uoAgAj9GgCSQL8ujkxh9SsLAL2gqqoqZs+eHcOGDYuIiMceeyyOPPLIaGhoiL59+8brr78e9913X0ycOLHElQJA5dKvAaD86dfFZaoWAHrVscceG/X19e3XJ0yYEK+99lr89Kc/jYULF8ZBBx0UO++8cwkrBAD0awAof/p1cRlxDgAAAAAAHZjjHAAAAAAAOhCcAwAAAABAB4JzAAAAAADoQHAOAAAAAAAdCM4BAAAAAKADwTkAAAAAAHQgOAcAAAAAgA4E5wAAAAAA0IHgHAAAAAAAOhCcAwAAAABAB4JzAAAAAADo4P8HAN1pHSVPps0AAAAASUVORK5CYII="
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "methods_pruned = [method for method in methods if method != 'HF']\n",
    "bases_pruned = ['aug-cc-pVTZ', 'aug-cc-pVQZ', 'aug-cc-pV[34]Z']\n",
    "fig, panes = plt.subplots(nrows=1, ncols=len(bases_pruned), sharey=True, figsize=(18, 6))\n",
    "\n",
    "for pane in range(len(bases_pruned)):\n",
    "    data = []\n",
    "    for method in methods_pruned:\n",
    "        data.append(\n",
    "            pymolpro.database.analyse(results[method][bases_pruned[pane]],\n",
    "                                      db,'kJ/mol')['reaction energy deviations'].to_numpy()[:, 0]\n",
    "        )\n",
    "    panes[pane].violinplot(data, showmeans=True, showextrema=True, vert=True, bw_method='silverman')\n",
    "    panes[pane].set_xticks(range(1, len(methods_pruned) + 1), labels=methods_pruned, rotation=-90)\n",
    "    panes[pane].set_title(bases_pruned[pane])\n",
    "panes[0].set_ylabel('Error / kJ/mol')\n",
    "plt.savefig(project_name + \".violin.pdf\")"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-02T06:43:57.900086Z",
     "end_time": "2023-05-02T06:43:58.105675Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/var/folders/qd/f26v1rbd3g5341rsg3s9j4v80000gn/T/ipykernel_71408/3833976484.py:1: FutureWarning: In future versions `DataFrame.to_latex` is expected to utilise the base implementation of `Styler.to_latex` for formatting and rendering. The arguments signature may therefore change. It is recommended instead to use `DataFrame.style.to_latex` which also contains additional functionality.\n",
      "  basis_errors.to_latex(float_format='%.2f')\n"
     ]
    },
    {
     "data": {
      "text/plain": "'\\\\begin{tabular}{lrrrrr}\\n\\\\toprule\\n{} & \\\\multicolumn{5}{l}{DF-LMP2} \\\\\\\\\\n{} & aug-cc-pVDZ & aug-cc-pVTZ & aug-cc-pVQZ & aug-cc-pV[23]Z & aug-cc-pV[34]Z \\\\\\\\\\n\\\\midrule\\nMSD    &       -1.21 &       -4.05 &       -0.54 &          -4.85 &           2.07 \\\\\\\\\\nSTDEVD &        1.36 &        5.71 &        0.26 &           7.71 &           3.80 \\\\\\\\\\nMAD    &        1.44 &        4.05 &        0.54 &           4.85 &           2.24 \\\\\\\\\\nMAXD   &        2.23 &       10.62 &        0.79 &          13.76 &           6.46 \\\\\\\\\\nRMSD   &        1.64 &        6.17 &        0.58 &           7.95 &           3.73 \\\\\\\\\\n\\\\bottomrule\\n\\\\end{tabular}\\n'"
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "basis_errors.to_latex(float_format='%.2f')"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-02T07:01:12.843041Z",
     "end_time": "2023-05-02T07:01:12.847693Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\\begin{tabular}{lrrrrr}\n",
      "\\toprule\n",
      " & \\multicolumn{5}{c}{DF-LMP2} \\\\\n",
      " & aug-cc-pVDZ & aug-cc-pVTZ & aug-cc-pVQZ & aug-cc-pV[23]Z & aug-cc-pV[34]Z \\\\\n",
      "\\midrule\n",
      "MSD & -1.210 & -4.046 & -0.539 & -4.854 & 2.070 \\\\\n",
      "STDEVD & 1.362 & 5.710 & 0.262 & 7.711 & 3.804 \\\\\n",
      "MAD & 1.435 & 4.046 & 0.539 & 4.854 & 2.238 \\\\\n",
      "MAXD & 2.226 & 10.624 & 0.792 & 13.757 & 6.462 \\\\\n",
      "RMSD & 1.644 & 6.173 & 0.580 & 7.950 & 3.733 \\\\\n",
      "\\bottomrule\n",
      "\\end{tabular}\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print(basis_errors.style.format(precision=3).to_latex(hrules=True,multicol_align='c'))"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-05-02T07:08:32.681098Z",
     "end_time": "2023-05-02T07:08:32.683913Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "outputs": [],
   "source": [
    "with open(project_name + '.tex', 'w') as tf:\n",
    "    tf.write('\\\\def\\\\toprule{\\\\hline\\\\hline}\\n\\\\def\\\\midrule{\\\\hline}\\n\\\\def\\\\bottomrule{\\\\hline\\\\hline}')\n",
    "    tf.write(df_exp_reaction_meanerror.style.to_latex())\n",
    "    tf.write(df_exp_reaction_std.style.to_latex())\n",
    "    tf.write(df_exp_reaction_meanabserror.style.to_latex())\n",
    "    tf.write(df_exp_reaction_maxerror.style.to_latex())"
   ],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-04-12T09:44:10.342725Z",
     "end_time": "2023-04-12T09:44:10.346654Z"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "outputs": [],
   "source": [],
   "metadata": {
    "collapsed": false,
    "ExecuteTime": {
     "start_time": "2023-04-12T09:44:10.345111Z",
     "end_time": "2023-04-12T09:44:10.346749Z"
    }
   }
  }
 ],
 "metadata": {
  "kernelspec": {
   "name": "python3",
   "language": "python",
   "display_name": "Python 3 (ipykernel)"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.6"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}